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We apply the effective potential analytic continuation (EPAC) method to the calculation of real 



in 
o 

time quantum correlation functions involving operators nonlinear in the position operator q. For 
a harmonic system the EPAC method provides the exact correlation function at all temperature 
^-j, ranges, while the other quantum dynamics methods, the centroid molecular dynamics and the ring 

polymer molecular dynamics, become worse at lower temperature. For an asymmetric anharmonic 
system, the EPAC correlation function is in very good agreement with the exact one at t = 0. 
When the time increases from zero, the EPAC method gives good coincidence with the exact result 
at lower temperature. Finally, we propose a simplified version of the EPAC method to reduce the 
■ computational cost required for the calculation of the standard effective potential. 
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O 
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in 
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The imaginary time path integral [1J has provided a useful framework suitable for numerical analyses of quantum 
O |' statistical-mechanical systems. Most of the static properties of quantum systems can be calculated by means of the 
path integral Monte Carlo (PIMC) or path integral molecular dynamics (PIMD) technique 0, However, it is 
not straightforward to apply the PIMC/PIMD methods to computing dynamical properties such as the real time 
quantum correlation function (A(t)B(0))p. This is because it is nontrivial to construct real time quantities from a 
finite number of imaginary time data obtained numerically |4j . To overcome such difficulty, a number of promising 
methods of numerical analytic continuation based on the maximum entropy method have been proposed and applied 



to various many-body systems 
1%, 





Recently, the centroid molecular dynamics (CMD) method 9], the ring polymer molecular dynamics (RPMD) 
method jlCj . and the effective potential analytic continuation (EPAC) method 'llj have been proposed as new quantum 
dynamics methods to calculate real time quantum correlation functions at finite temperature. Both the CMD and 
the RPMD are the methods to calculate the canonical (Kubo-transformed) correlation function 01 (A(t) B (O))^* 1 by 
means of molecular dynamics techniques. On the other hand, the EPAC is a method to obtain the real time quantum 
correlation function (A(t)B{Q)) p by means of the effective action formalism and an analytic continuation 

procedure 0| . It has been shown analytically that all these methods are exact in harmonic systems for the real time 
quantum correlation functions of a linear function of the position operator, A = B = q |lfl Il7l Il8|. However, for 
nonlinear operators such as A = B = q n , it is nontrivial whether these quantum dynamics methods yield the exact 
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result even in harmonic systems . This is the nonlinear operator problem in quantum dynamics methods. 

From a practical point of view, the nonlinear operator problem is a quite important subject to be tackled. In 
many problems of chemical interest, the real time correlation functions of nonlinear operators are often needed for 
the calculation of various dynamical properties 0]. For example, the CMD has been applied to the calculation of the 
transport coefficients such as thermal conductivity, shear viscosity, and bulk viscosity of quantum liquid parahydrogen 
|'2p| . Here it is found that the calculated transport properties are in good agreement with the experimental data. 
However, there is no rigorous theoretical basis for applying the CMD method to such properties represented by means 
of the correlation functions involving operators nonlinear in q (or momentum operator p). Therefore we need, in 
general, a quantum dynamics method which is theoretically valid even for the time correlation functions of nonlinear 
operators. 

The present status of the nonlinear operator problem in the three quantum dynamics methods, CMD, RPMD, 
and EPAC, is summarized as follows. For the CMD method, Reichman et at have argued this problem in their 
pioneering paper to conclude that a CMD correlation function involving nonlinear operators corresponds to a higher- 
order Kubo-transformed correlation function 21] . On the other hand, Craig and Manolopoulos have shown that the 
RPMD is exact for all the operators involving q in the limit t — > 0]. As for our EPAC method, the nonlinear 
operator problem has not been examined yet. In addition to these three methods, a theoretical approach based on 
the quantum mode-coupling theory has been applied to study the dynamical properties involving nonlinear operators 
in quantum liquids [sib^. 

In the present paper, we develop the method of the EPAC for the calculation of the real time quantum correlation 
function involving the nonlinear operator A = B = q n . As a simple example, at first we show the EPAC correlation 
function (q 2 (t)q 2 (0))^ PAC for a harmonic oscillator comparing with the results of the CMD and the RPMD. Next we 
calculate (q 2 {t)q 2 (0))p PAC numerically in an asymmetric anharmonic system. We also propose a simplified EPAC 
method to reduce the computational cost required in the EPAC calculation. 

In Sec. II, we summarize the effective action formalism and present how to calculate the EPAC correlation function 
involving nonlinear operators. The results for a harmonic oscillator are given in Sec. III. Numerical results for 
an anharmonic oscillator are shown in Sec. IV. In this section we also present the simplified EPAC method. The 
conclusions are given in Sec. V. 



II. EFFECTIVE POTENTIAL ANALYTIC CONTINUATION METHOD FOR NONLINEAR 

OPERATORS 



Hereafter we treat the real time quantum autocorrelation function of the nonlinear operator q n , (q n (t)q n (0))p. 
In principle, this can be obtained from the imaginary time Green function (Tq n (r)g™(0))^ via an analytic continuation 
procedure [lfij . Here T represents a time-ordered product. On the other hand, it is known that the Green function 
(Tq n (T)q n (0))i3 is given as a special case of 2n-point imaginary time Green function {Tq(ji) ■ ■ ■ q(T?.r,.))iH : which can 
be constructed from the standard effective potential Vg(Q) appearing in the effective action formalism \lA ll-d Il5|. 
Consequently, the real time quantum correlation function (q n (t)q n (0))p should be obtained by means of the effective 
action formalism and the analytic continuation. A series of these procedures is the EPAC method for the nonlinear 
operator q n , which we newly show in this section. As the simplest example, we present the EPAC calculation of the 
autocorrelation function of the quadratic operator q 2 , {q 2 (t)q 2 (0))/3, in Sees. IIA and IIB. 
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A. Effective action formalism for multipoint Green functions 



We begin with the effective action formalism |13t Il4l Il5| . Consider a quantum system where a quantum particle 
of mass m moves in a one-dimensional potential V(q) at inverse temperature (3. The quantum canonical partition 
function of this system is expressed in terms of the imaginary time path integral 



Zp[J] 



rq(f3h)=<l 

dq / T>q exp 

Jq(0)=q 



-s E in- 



dTJ{r)q(r)lh 



mq 2 + V(q) 



where Se is the Euclidean action functional 

,/jft 

Sslq] = I dr 
Jo 

and J(r) is an external source. The generating functional W^[J] is defined as 

Wp[J] = hlogZp[J]. 



(1) 



(2) 



(3) 



The functional derivative of W^[J] with respect to J(t) produces the quantum statistical-mechanical expectation 
value of the operator q in the presence of J, 



7j(7T = {q{T)) ? = 0(r) - 



The effective action Tp[Q] is defined by the Legendre transform of W^J], 

r (3h 



■Wp[J] 



drJ(T)Q{r), 



(4) 



(5) 



which satisfies the quantum-mechanical Euler-Lagrange equation 5Tp[Q]/SQ(r) = J(r). The exact quantum 
statistical-mechanical expectation value Qs(t)(= (<?(t)}^) is obtained as a solution of the equation 



jTgjg] 
5Q(t) 



(6) 



Note that the expectation value Qs(t) is independent of imaginary time r in thermal equilibrium. 

The n-point connected Green function G^(ti, ...,t„) is generated by the nth-order functional derivative of W^g[J] 
with respect to J(r) 0,0]. For example, the two-point connected Green function is given by 



5J(ti)5J{t 2 ) 



J=o 



= {Tq{Ti)q{T 2 ))p - (q(n)) (q(r 2 ))p. 

Then we obtain the two-point Green function 

(Tq( Tl )q(T 2 ))p = G%(n - r 2 ) + Q S (ti)Q S (t 2 ). 
The procedure in Eqs. Q-© has been described in Ref. 0. In a similar way, the three- and four-point Gr 



(7) 



(8) 



ecn 



4 



functions can be explicitly expressed as 



7=0 



(Tq(Ti)q(r2)q(T 3 )q(r 4 ))p 



SJ( Ti )6J(t 2 )6J(t 3 ) 
+Qs(Ti)(Tq(T 2 )q(T 3 ))f 3 + (1 

- 2 Qs(ti)Qs(t 2 )Q s (t 3 ), 



2) + (!<-» 3) 



5J(ti)SJ(t 2 )SJ{ts)SJ(t4) 
S 3 W fJ [J] 



,7=0 



7=0 



+H 0s(Tl) ^J(r 2 )5J(r 3 )5J(r 4 
+ (Tq(T 1 )q(T 2 )) p (Tq(T 3 )q(r 4 )) p + (2 

- 2 Qs(ti)Qs(t2)Qs(t 3 )Qs(t4), 



(1-2) 
3) + (2 



(1 <-► 3) + (1 4) 
>4) 



(9) 



(10) 



respectively. Here (i j) denotes the interchange of Tj with Tj. 

On the other hand, the nth-order functional derivative 5 n Wp/5J n is connected with the nth-order functional 
derivative of the effective action 8 n T^/SQ n 



S 2 Wp[J] 



5. ] (n)S. ] (r 2 ) 



<W(ti)<5J(t 2 )<5J(t 3 ) 



,7=0 



.7=0 



<JQ(ri)<SQ(r 2 ) 



-«s / ^ 



/3fi 

ds 2 



/3ft 

ds 3 



xG^(si - ri)G^(s 2 - T 2 )G ( f ) {s 3 ~ t 3 ) 

8 3 rp[Q] 



5Q( Sl )5Q(s 2 )6Q(s 3 ) 



S 4 W,[J] 



5J{ri)5J{T2)8J{r 3 )5J{Ti) 



.7=0 



8h 
ds 2 



/3ft 

ds 3 



pn 

d,S4 



xG%( Sl - t 1 )G c (s 2 - T 2 )G c p {s 3 - t 3 )GC(s 4 - r 4 ) 



x 



6Q( Sl )5Q( S2 )6Q(s 3 )6Q(s 4 ) 



1 



Ph 

dsi 







II 



/3ft 

ds 2 



/3ft 

ds 3 



ph 
ds 4 



pn 
du 



8h 
dv 



xCp ( ai - T1 )G% (s 2 - r 2 )G^ ( M - t;)G^ (s 3 - t 3 )G£ (a 4 - r 4 ) 



6Q( Sl )5Q(s 2 )6Q(u) 
(2 «-> 3) + (2 <-> 4). 



<5 3 I>[Q] 



q=q s <JQ(u)5Q(s 3 )(SQ(a4) 



(11) 



(12) 



(13) 



Using Eqs. (|7I)- (|13|I . we can construct the imaginary time Green functions in terms of the effective action r^[Q] 
Now we employ the local potential approximation (LPA) to the effective action, 



/3ft 



dr 



Vp(Q) + ^mQ 2 



(14) 



where Vp(Q) is the standard effective potential, i.e., the leading order of the derivative expansion of [T3L I14I IT^ . 

For the case Qs{t) is r independent, Eq. JJjJ becomes the stationary condition ^^/^ < 3|g = Q min = 0> which determines 
the expectation value Qs{t) as the standard effective potential minimum Qmin- Then the functional derivatives of 
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the effective action 5 n Tp/SQ n become 



■Pi 



SQ( Tl )SQ(T2) 
S 3 Tp[Q] 



d 2 v p 



dr 2 dQ 



?=Q„ 



<5Q(ti)<5Q(t 2 )5Q(t 3 ) 
S 4 Tp[Q] 



d 3 Va 



dQ 3 
d 4 Vp 



SQ{tx)SQ{t2)8Q{t 3 )8Q{t a 
If we expand Vp around the minimum Q = Q 



dQ 4 



S(n - r 2 ), 
S(n - t 2 )S(ti - r 3 ), 
S(n - t 2 )S(t 1 - t 3 )S(ti - r 4 ). 



W) = E^7(3-QW> n - 

i n! 

n=0 

then the derivatives of appearing in Eqs. (|15(1 - (|17[1 are given as the coefficients in the series, 
d n Vp/dQ n \ Q=Qmin . 



(15) 
(16) 
(17) 

(18) 



B. EPAC method 



Now we proceed to the calculations of the real time quantum correlation function (q 2 (t)q 2 (0))p, i.e., the aim of 
the present paper. First, by setting t\ = t 2 — r (0 < r < (3h) and T3 = T4 = in Eq. I|l()|l. we obtain the imaginary 
time Green function 



(Tq 2 (r)q 2 {0))p 



S 4 Wp[J] 



5J(t)8J{t)6J(0)6J(0) 



j=o 



2h 2 Q z 



8 3 W P [J] 



5J{t)8J(t)5J(0) 



+ 



5 3 Wp[J] 



=0 6J(t)6J(0)6J(0) 



j=o. 



2 ((Tq( T )q(0))f>) 2 + ((T<j(0)g(0))^) 2 - 2 Q 



4 

min ' 



(19) 



With the definition ujp = y/ 1 a 2 jm and a — (3fkvp/2, from Eqs. Q, ©, l{TT jl -l{T3 )l . and l(T^ )l -l(T8 ^ . the components of 
Eq. (|19(l can approximately be written as 



(Tq( T )q(0))p 



2mu!p e a — e 



) 2 . 

' mm 5 



(20) 



S 3 Wp[J] 



5J(t)5J(t)6J(0) 



6 3 Wp[J] 



J=0 SJ(t)5J(0)5J(0) 
a 3 1 



,7=0 



6m 3 a^ ( e Q - e - Q ) E 
+ 2(e 



[(e Q 



3a _ -a\ „-WpT 



+ 2(e a 



3a\ „ujqt 



-3a _ e -a^ e 2w^T 

-6(e Q -e-«)] , 



(21) 



5 A W fj [J] 
5J(t)6J(t)6J(0)5J(0) 



j=o 

-^-g j [ 4a(4 + e~ 2 ^ + e 2 ^ ■ 



32m 4 c^ ( e a _ e -a) 4 
+2(e 4a - l^re" 21 ^ + 2(e~ 4Q - l)w /3 re 2 ^ r 
+8(e 2a - e~ 2Q ) + (e 4Q - l)e- 2 ^ T + (1 - e- 4Q )e 2 ^ T 

n 2 1 



288m 5 wl ( e o — e~ Q ) 



[ 60a(e Q - e^ Q )(4 + e~^ T + e 2 ^ T ) 



+30(e 5Q - 


e 3a - 


e Q +e~ Q V /3 re- 2 ^ r + 30(- 


-e Q + e~ a + e- 3a - e'^ujpTe 2 "^ 


+16(e 5Q - 


2e 3a - 


fe a )e _ ^ T + 16(e" a - 2e' 


3a _|_ e -5a^ e 3aj/3T 


-17(e 5Q - 


e 3a - 


e a + e- a )e~ 2u} ^ - 17(e Q - 


er a - e~ 3a + e - 5Q )e 2 "" T 


+16(e 5Q + 


3e 3Q - 


- 8e Q + 3e~ Q + e~ 3Q )e^ T 


+ 16(e 3a + 3e Q - 8e" Q + 3e~ 3Q + e 


+248(e 3Q 


-e Q - 


e - Q+e -3a- ) ] _ 





~ 5Q )e"' 3T 



(22) 



Substituting Eqs. (|20Jl - 122|l into Eq. 1)19(1 . we can approximately represent the quadratic imaginary time Green 
function (Tq 2 (T)q 2 (0)) p in terms of the quantities Q mul , up, as, and 04. 

Then, by means of the analytic continuation the EPAC correlation function (q 2 (t)q 2 (0))® PAC is obtained from 
the imaginary time quantity (Tq 2 (T)q 2 (0)} p. This is the final step of the EPAC method for the correlation function 
of the nonlinear operator q 2 . The EPAC correlation function, which is the approximation to the exact correlation 
function {q 2 {t)q 2 (O)) , can be expressed in a simple form 



where 



(q 2 (t)q 2 (0))^ = a 4 A(t) + a 2 3 B(t) + a 3 Q min C(t) + D{t), (23) 



Ait) = - Z £ H {ea _\_ a)i [ 8(e 2 " - e- 2 ") + 8a(cos2^ + 2) 
+ (e 4a - e~ 4a ) {2u)pt sin 2u>pt + cos 2u)pt) 

+ i (e 4a + e~ ia - 2)(2uj tco8 2(jf3t - sin 2^) ] , (24) 

B{t) = 28&4 (e° =W [ 248(e3 " " 6 " ~ 6 "° + e ~ 30) + 12 ° Q(e " " e '")( COs2 ^ t + 2 ) 
+16(e 5Q - 2e 3Q + e a + e~ Q - 2e~ 3 " + e~ 5Q ) cos 3^* 
+ (e 5Q - e 3a - e~ 3a + e~ ba )(iQujptsai2upt ~ 17 cos 2^3*) 
+16(e 5Q + 4e 3a - 5e a - 5e~ Q + 4e~ 3a + e" 5 ") cos^i 
+ i {-16(e 5Q - 2e 3Q + e a ~ e~ a + 2e~ 3Q - e~ 5Q ) sin 3^ 

+ (e 5Q - e 3a - 2e a + 2e~ a + e~ 3Q - e' 5a )(30ujf 3 tcos2upt + 17 sin 2^) 

-16(e 5Q + 2e 3Q - lie" + lie"" - 2e~ 3a - e~ 5a ) sinwpt } ] , (25) 

C (t) = " Tr[6(e Q -e- Q ) 

3m 3 w^ (e Q - e~ Q ) 3 L 

-(e 3Q - e Q + e" Q - e" 3 ") cos 2cu p t + 2(e 3a + e a - e~ Q - e" 3 ") cosupt 

+ i (e 3a - e a - e' a + e" 3Q )(sin2w (3 t -2saxw p t) ] , (26) 

h 2 

D(t) = 5 — 5- T 2cotha(coth2acos2o> i at — isxn2u>gt) + 2coth a— ll 

+ min [ 2cothacosw /3 i - i 2 sin o;^ + coth a } + Qf nin . (27) 
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It should be noted that this correlation function, Eq. (|23|) . consists of three oscillation modes with the frequencies 
u>p, 2u>p, and Stop; the coefficients 03 and 04 contribute only to the amplitude of the oscillations. 

Similarly, the higher-order real time correlation functions (q n (t)q n (0))p can approximately be obtained as 
(<?"(£) 9" (0)}Jp AC using the information of the standard effective potential Vp{Q). 



III. SIMPLE EXAMPLE : A HARMONIC SYSTEM 



In this section, as a simple example we show the calculation of the EPAC correlation function {q 2 (t)q 2 (0))Jp AC for 
a quantum harmonic oscillator whose classical potential is given by 

V{q)= l -mu?q 2 . (28) 



Next, the results of the other quantum dynamics methods, the CMD and the RPMD [ljj, are shown 



EPAC correlation function for a harmonic oscillator 



The standard effective potential for the system l|28[) is obtained as [l 

V (0) = \mu?Q 2 + i log ^2 sinh ^\ . (29) 

For the harmonic system, the minimum of Vg(Q) is located at the point Q = Q mm = 0, while the effective frequency 
is top — lu. Furthermore, it is evident that the higher-order derivatives of Vp(Q) vanish, 03 = 04 = 0. Consequently, 
the EPAC correlation function H23[) becomes 



EPAC 



Am 2 u 



2 coth — — (coth [3huj cos 2ujt — i sin 2u>t) + 2 coth 2 — 1 



(30) 



which is equal to the exact quantum correlation function (g 2 (t)q 2 (0))^. That is, the EPAC method is exact in the 
harmonic system l|28|l at any temperature. This is because the LPA [Eq. I|14jl]. which is the only approximation 
employed in the EPAC method, is exact for harmonic systems [2^. Therefore, it can be shown that the EPAC 
correlation function (q n (t)q n (0))® PAC is exact for any n for harmonic systems. 



B. The other quantum dynamics methods 



In harmonic systems, both the CMD and the RPMD are exact for the linear operator 0,0]. That is, both the 
CMD correlation function and the RPMD correlation function are equal to the exact canonical, or Kubo-transformed, 
correlation function {q(t)q(O)} c an = (1//3) Jjf dX (q(t-ih\)q(O)) . However, even in harmonic systems, these methods 
are no longer exact for the correlation functions of nonlinear operators Il7| . 

In the harmonic system l|28ll . we explicitly present the exact canonical autocorrelation function of the nonlinear 
operator q 2 



(q 2 {t)q 2 mT = 



■ coth cos 2ujt + 2 coth 



Jhu> 2 2 

while the corresponding CMD correlation function with the classical operator Q is obtained as 



(31) 



(q 2 c (t)q 2 c(0))p MB = -^^(0082^ + 2), (32) 



where q c is the position centroid variable. The CMD correlation function with the effective classical operator |l7| is 
also obtained as 



(&w)m)p = ^Kj~^ cos2Lut + — coth — + ^ ^ 

where {q 2 )p is the effective classical operator: (q 2 )p = q\ + (7i/2mw)(coth(/37Kj/2) — 2/(3fkS). On the other hand, the 
RPMD correlation function is given by 



«W(0)>r MD = Mm 



fPm 2 



P 1 P P 1 

V -r (cos 2w„i + 1) + V V -5-3 

n=l ™ n=l Z=l 71 ' 



(34) 



where uo n = \/lu 2 + (2fcp/m)(l — cos(27rn/P)), fcp = mP 2 /(f3 2 h 2 ), and P is the number of discretization of imaginary 
time. 

Figure ^ shows the plot of the four correlation functions, Eqs. <|31|) - (|34ll . at two different temperatures with the 
parameters h=kB = i rn = uj = l. Here we have also set P as 1000, which is so large as to make the RPMD results 
converged sufficiently. At higher temperature (3 = 1, in Fig. n( a )i au the CMD results and the RPMD result are 
slightly deviate from the exact one. These deviations become remarkable at lower temperature [3 = 10 [Fig. [2(b)]. 
As is observed clearly at (3 = 10, the CMD correlation functions [Eqs. (|32|l and l|33|) ] fail to reproduce the exact value 
at t = and the exact amplitude, while they oscillate with the correct frequency 2ui. On the other hand, the RPMD 
correlation function Ij34(l is exact at t — 0, while it damps with time because of the dephasing effect caused by the 
mode summation y~) j . This damping is remarkable at lower temperature. 

For the CMD method, Jang and Voth have explicitly shown that the CMD should be used for the computation 
of the canonical correlation function (A(t)B(0) )T^ = {1//3) d\ (A(t — ih\)B(0))p where the operator A must be 
linear in position and/or momentum operators |17| . As for the operator B, the effective classical operator B^ should 
be used for the CMD description of this correlation function 17]. The CMD method is known to be problematic 
when the operator A is a nonlinear function. Afterwards two major approaches have been reported to deal with this 
nonlinear operator problem in the CMD method [21I l^. Here we note these approaches: 

(1) The first approach is a theory based on the higher-order Kubo-transformed correlation function; Reichman et 
al. have pointed out that the CMD correlation function with the effective classical operator Bp, (q™(i)i?|(0)}^ MD , 
corresponds to the nth-order Kubo-transformed correlation function (q n (t) B (0))^ th j^]. That is, the CMD correlation 
function with effective classical operator, {q 2 {t){q 2 ) c p{Q)) ( ^ WD , corresponds to the second-order Kubo-transformed 
correlation function 



2 



Jo 



(r(t)q 2 (0))p nd = I d\ I dq (q(t-ihX)q(t-ihr,)q 2 (0))p. (35) 



For the harmonic oscillator (|28|l . it can readily be shown that Eq. 133|) is equal to Eq. (|35|l . Thus, the CMD 
with effective classical operators works very well in this framework. However, there remains a practical problem in 
converting the higher-order Kubo-transformed correlation function (q n (t)B(0)y^ th to the original quantum correlation 
function (q n {t)B(Q))p; this conversion would be a complicated procedure in general. 0| 

(2) The second approach is based on novel expressions for physical quantities, which are formulated using correlation 
functions involving linear operators. Geva et al. have shown that the quantum reaction rate constant can be expressed 
in terms of the canonical correlation function (q(t)B(0)) o an |2J|. Their approach enables us to evaluate the rate 

nTHr 



constants via CMD calculations without further approximations [2J, |25| 



9 



TABLE I: The results of the quantities for the quantum anharmonic oscillator 136> , 



p 




U)0 


a-i 


(J4 


0.1 


-1.3735019 


1.07083695 


0.10132291 


0.1018375 


1 


-0.3375973 


0.91069063 


0.41549732 


0.3305302 


10 


-0.1501482 


0.96628105 


0.54407872 


0.2606658 


100 


-0.1501276 


0.96631313 


0.54396628 


0.2608735 



IV. NUMERICAL TESTS : AN ANHARMONIC SYSTEM 



In this section we calculate the EPAC correlation function (q 2 (t)q 2 (0))^j PAC for an asymmetric anharmonic 



system with the classical potential 0, 0, 0, Q, 0] 

V{q) = \q 2 



10 



100 



(36) 



where natural units h = ks = m = 1 are employed. In the following, first we show the calculated EPAC correlation 
functions with the exact correlation functions (q 2 {t)q 2 (0))p, and then we present a simplified version of the EPAC 
method to discuss its validity. 



A. EPAC results for an anharmonic oscillator 

To evaluate the standard effective potential Vp{Q), we need to compute the generating functional W/3[J] [Eq. 
0] and to carry out the Legendre transformation (JSJ). There are various computational schemes to evaluate the 
generating functional Wp [J] . Among them, the PIMD /PIMC technique is a practical tool to evaluate Wp [J] directly 
or indirectly [lli lisi . Here we have, however, employed the renormalization group (RG) method, which is suitable 
for precise calculation of VF|a[J] j^]. Then the standard effective potential Vp(Q) has been computed by means of the 
numerical Legendre transformation Figure [2] shows the evaluated standard effective potentials at various inverse 
temperatures (3 — 0.1, 1, 10, and 100. Minimizing the effective potential Vp(Q), we have obtained the quantities Q m in> 
ujp, 03, and 04. Table^lists the results of computed quantities for the system (|36f) . 

Then we have obtained the EPAC correlation function (<7 2 (i)<? 2 (0))Jf PAC using the quantities listed in Table [I] For 
reference, we have also calculated the exact quantum correlation function (q 2 (t)q 2 (0))p by solving the Schrodinger 
equation numerically [2^. Figure |21 shows the real part of {q 2 {t)q 2 (0))® PAC together with the real part of {q 2 {t)q 2 (0)} p 
at various inverse temperatures (3 — 0.1, 1, and 10. At t — 0, each EPAC correlation function is in very good agreement 
with the exact correlation function. This means that the LPA [Eq. H14fl ] , which is the only approximation employed 
in the EPAC method, is fairly good for the calculation of the static property (q 2 (0)q 2 (0))p. In fact, good results have 
also been obtained in the calculation of the EPAC correlation function involving the linear operator q, (<j(0)(j(0))Jf PAC , 
for the system ^ Q- 

On the other hand, as time t increases, each EPAC correlation function deviates from the exact one and this 
deviation becomes worse at higher temperature. There are a couple of reasons for such deviation. The first is the 
number of oscillation modes in the correlation functions. Although the EPAC correlation function H23|l consists of a 
limited number of oscillation modes with the frequencies up, 2w^, and 3w^, the exact correlation function (q 2 (t)q 2 (0)) p , 
in general, consists of many oscillation modes. Then at higher-temperature (3 = 0.1 [Fig. 0(a)], a larger number 
of oscillation modes get to contribute to the exact correlation function, resulting in rapid damping. The EPAC 
correlation function (q 2 {t)q 2 (0)}^ PAC including only three oscillation modes cannot represent such damping behavior. 
On the other hand, as the temperature lowers [Fig. [3] (b) and (c)], the EPAC correlation function becomes very 
closer to the exact one because oscillation modes involved in the exact correlation function become fewer. This is 
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the same behavior as we observed in the calculation of (q{t)q(0))® PAC for the system iKSOt [is| . The second reason 
for the disagreement between the EPAC and the exact results at long time is found in the anomalous terms in the 
EPAC correlation function (|23|) . In Eq. (|23p. the components A(t) and B(t) both contain the terms proportional 
to tsm2u>pt or tcos2copt, which exhibit amplified oscillations and diverge in the long time limit. It is expected that 
these anomalous terms disappear only if we calculate the higher order derivative terms in the effective action r^[Q] 
beyond the LPA [Eq. I)14|l] . A divergence- free EPAC correlation function could also be obtained if we omitted all the 
terms proportional to re -2 "* 31 " or Te 2w/3T in Eq. (|22(l . However, it should be noted that such prescription would break 
the periodicity of the imaginary time Green function, (Tq 2 ((3h)q 2 (0))fj — {Tq 2 (0)q 2 (Q)) p. 

Finally we mention how the EPAC method could capture the anharmonic effects in quantum statistical systems. In 
the EPAC method, all the quantum/thermal effects are included via the standard effective potential Vp(Q). For the 
harmonic case l|28|l . the effective frequency equals the classical frequency (a-> c ff = uj), and the higher-order coefficients 
vanish (a„ = for n > 3). However, for anharmonic systems, the effective frequency w e ff deviates from the classical 
one and the higher-order coefficients a„>3 have nonzero values. As seen in Eq. (|23|) . the coefficients a„>3 contribute 
only to the amplitude of the oscillations of the EPAC correlation functions, while the effective frequency w e ff provides 
oscillation modes with the frequency lujp (I: integer). The EPAC is a method to approximate the exact dynamics 
of quantum systems using a finite number of oscillation modes with the frequency lujp, and therefore the EPAC 
correlation function always exhibits harmoniclike oscillations. That is, although the anharmonic effects are included 
via the quantities w e g and a„>3, the EPAC correlation functions cannot reproduce a certain type of anharmonic 
effects such as the dephasing and the rapid damping at higher temperature [Fig. 01(a)]. This quasiharmonic property 
of the EPAC correlation function comes from the LPA [Eq. (|14fl ]. the only approximation employed in the EPAC 
method. At the same time, we should note an important problem arising in the standard effective potential approach 
when it is applied to quantum systems that include the dissociation limit lim^i^^ V(q) = |18l l29j: for example, 
the systems represented in terms of the Morse potential V(r) (r = \q± — q2\) are classified into this category. For such 
systems with the dissociation limit, the minimum should disappear in the standard effective potential Vg(i?) because 
it is always convex for < R < oo according to its definition [ll], |lj, |lj| . Then the frequency W/j cannot be defined; 
this would lead to some unphysical flaw. Therefore, the EPAC is not very suitable for the full description of, e.g., the 
dissociation reaction of isolated diatomic molecules. Still we note that at low temperatures the EPAC does work well 
because the po lyno mial expansion around the potential minimum r m ; n is so valid as to neglect the dissociation limit 
approximately 18]. 

As for quantum liquids whose intermolecular interaction is represented in terms of, e.g., the Lennard- Jones potential 
with the dissociation limit as well, the situation of such many-body systems is not so simple as an isolated diatomic 
system, because potential-energy surface in configurational space includes many local minima [i(3, 13 • Reichman 
and Voth have discussed various effective harmonic theories for liquid dynamics associated with the curvature of such 
potential-energy surface |Q| . Certainly, there is a similarity between the effective harmonic theories for liquids and the 
standard effective potential approach we are discussing, in that ujp is defined as the second derivative of the standard 
effective potential Vp(Q) as well. However, the standard effective potential and the derived EPAC method is not 
suitable for describing the diffusion in liquids because some well-defined multidimensional Vp(Q) should be convex to 
have a single minimum for one-dimension, resulting in just the oscillatory motions with real positive frequencies u>$ . 
Thus, at the present stage, the EPAC method does not fully capture the molecular diffusion in quantum liquids 33]. 

Rather, the EPAC method should be useful for the investigation of, e.g., the coherent dynamics in bound systems 
involving proton transfer reaction [34] ] for which the model potentials are typically represented as double-well type 
such as V(q) = —aq 2 + bq 4 (a, b > 0). These potentials are asymptotically superlinear, limigi^oo V(q)/\q\ > 0, without 
dissociation limit. In such potential systems the quasioscillating behavior of a correlation function is essential. It is 
true that this oscillation is a consequence of the quantum coherence between the localized states in the local minima 
of V(q). The EPAC method based on Vp(Q) is suitable for describing such oscillation [n|, because its frequencies lu^ 
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are evaluated from the convex effective potential Vp(Q), which includes the effects of quantum interference between 
the localized states 

B. Truncated EPAC method 

In this subsection we present a further approximation scheme useful for practical computation. Here we treat 
the truncated expansion of Vp around the minimum Q — Q m i n , 

2 

W) = E^rW-Q™ n ) n . (37) 

n=0 

This means that we approximate the standard effective potential Vp as an effective harmonic potential and omit all 
the higher order derivatives of Vp(Q), a n = d n Vp / <9Q n |Q=Q min (n > 3). Since the determination of a„>3 requires a 
heavy computation of Vp(Q) with high precision, the truncation such as Eq. (|3Tjl reduces the computational cost 
significantly. 

Using Eq. 1(3 7|) , we obtain the truncated EPAC correlation function 

(eit)q 2 mf PAC = 0®, (38) 

where D(t) is just the function given in Eq. (|27p. This is the completely harmonic version of the EPAC method. 
Figure 0] shows the real part of Eq. I(38|) together with the real part of the exact quantum correlation function for 
the quantum anharmonic oscillator 1)36(1 . The temperatures are the same as in Fig. Here we have used again the 
results of Q m i n and ujp listed in Tabic [I] 

We see that, at t = 0, the truncated EPAC correlation functions are in good agreement with the exact results 
except for the middle temperature /? = 1 [Fig. |H(b)]. This behavior can be explained qualitatively as follows. Since 
cotha > 1, the highest power of cotha is regarded as the dominating factor in each term in the full EPAC correlation 
function J2S1) at t = 0, 



a 4 A(Q) -i-^fcoth 3 a, (39) 

SI Up 



a 3 fl (0) ~ 44 coth 'V (40) 
yo uj p 



a 3 Q min C(0) ~ -i^%icoth 2 a, (41) 

3 _L_ 

4w? 



D(0) ~ -— coth 2 a, (42) 



where Up ~ 0(1) (see Table [|J. For lower temperature (3 — 10, since cotha ~ 1 and 03,04 < 1 (see Table HJ, all 
the anharmonic terms in Eq. i)L'MI) at t = 0, i.e., 04^(0), a§i?(0), and a-iQ min C{Q), become negligible. Therefore 
the truncated EPAC result at t = 0, (<f (0)<f (0)}J, EPAC [Eq. JjgJ], is a good approximation to the full EPAC result. 
On the other hand, as the temperature increases, cotha becomes larger than unity, and then the anharmonic terms 
a4A(0) and a§B(0) are no longer negligible in the full EPAC correlation function because they contain coth 3 a 
and coth 5 a, respectively. Thus the truncated EPAC correlation function without the anharmonic terms fails to 
reproduce the full EPAC result at the middle temperature (3 = 1 [Fig. 0] (b)]. However, for higher-temperature 
j3 = 0.1, the coefficients 03 and 04 become very small and therefore the standard effective potential Vp(Q) has an 
effectively harmonic shape. Consequently, the truncation such as Eq. 1)37(1 becomes a very good approximation to the 
full potential ((18(1 . and therefore the truncated EPAC can reproduce the full EPAC result very well again. 
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As for the long time behavior, a discussion can be made as follows. The truncated EPAC correlation function 
(g 2 (t)^ 2 (0))^ p consists of two oscillation modes with the frequencies ujp and 2u>p, and it lacks the oscillation mode 
with frequency "Stop which exists in the full EPAC correlation function (q 2 (t)q 2 (0))^ PAC . Nevertheless, the truncated 
EPAC correlation function correctly reproduces the oscillation appearing in the full EPAC correlation function at 
each temperature (Figs. |21and2J). This is because the oscillation mode with frequency 2>uip exists only in B{t) and its 
contribution to Eq. \2'A\ is relatively small. It should also be noted that the truncated EPAC correlation function, only 
D(t), is free from an amplified oscillation expected in Eq. Il2:>ll and it never diverges in the t — > oo limit. We therefore 
expect that the truncated EPAC method, Eq. 1)38(1 . works well in low-temperature systems. Thus the truncated 
EPAC method can be useful for the practical calculation of the nonlinear dynamical properties of low-temperature 
condensed phase systems because it should reduce the computational cost. 

V. CONCLUDING REMARKS 

In this paper, we have focused on the nonlinear operator problem in quantum dynamics methods. At first, we 
have shown how to apply the EPAC method to the calculation of the correlation function (q n (t)q n (0))p and have 
given the EPAC correlation function (q 2 (t)q 2 (0))® PAC as an example. It has been shown that the EPAC method is 
exact in a harmonic system, while the other quantum dynamics methods, the CMD and the RPMD, fail to reproduce 
the exact results even in the harmonic system. 

Then we have applied the EPAC method to the asymmetric anharmonic system. We have seen that the EPAC 
correlation function at t = agrees well with the exact correlation function (<7 2 (0)g 2 (0)),g. As for the long time 
behavior, the EPAC result becomes better at lower temperature. These good properties suggest that the EPAC 
method can be a useful quantum dynamics method for nonlinear correlation functions of low temperature systems. 
We have also seen that the EPAC correlation function contains the terms which cause amplified oscillation. We 
suppose that these anomalous terms would disappear only if we improved the approximation beyond the LPA. 

Finally we have discussed the truncated EPAC method which is given by the truncation of the standard effective 
potential Vp(Q). We have tested this method in the same anharmonic system to which the full EPAC method has 
also been applied. Although the truncated EPAC method is not so good as the full EPAC method, it works fine at 
lower temperature. Since the truncated EPAC needs only the information up to the second order derivative of the 
standard effective potential Vp(Q), it should significantly reduce the computational cost. Therefore, from a practical 
aspect, we expect that the truncated EPAC is more suitable than the full EPAC for complex many-body systems at 
lower temperature. 

In the present work, the discussions have been restricted in the evaluation of the correlation function (q n (t)q n (0)){3. 
However, the multitime quantum correlation functions, e. g., ( q(t-\)q(t2)Q(ta)q(0))n, are also important property espe- 
cially in the context of the nonlinear optical spectroscopy [35| . It is a challenging task to apply the EPAC method to 
such calculation of multitime quantum correlation functions in near future. 
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FIG. 1: The plot of the exact canonical correlation function Eq. 131L the CMD correlation function with the classical operator 
(CMD-co) Eq. the CMD correlation function with the effective classical operator (CMD-eco) Eq. Ip^jl . and the RPMD 

correlation function Eq. 134H for the quantum harmonic oscillator 1281 . (a) at /3 = 1 and (b) at /3 = 10. 
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FIG. 2: The inverse temperature /3-dependence of the standard effective potential Vs(Q) for the asymmetric anharmonic 
classical potential 13611 . In this plot we set V/3(Qmin) = 0. 
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FIG. 3: The plot of the real part of the exact quantum correlation function {q 2 (t)q 2 (Q))fj and the real part of the EPAC 
correlation function (<? 2 (i)? 2 (0))^ PAC for the quantum anharmonic oscillator 136H . (a) at f3 = 0.1, (b) at = 1, and (c) at 
P = 10. 
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FIG. 4: The plot of the real part of the exact quantum correlation function {q 2 (t)q 2 (0)}p and the real part of the truncated 
EPAC correlation function (<7 2 (t)g 2 (0))| EPAC for the quantum anharmonic oscillator (a) at (3 = 0.1, (b) at (3 = 1, and (c) 

at /3 = 10. 



